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An optimized dose delivery system using computer graphics simula, 
tion techniques and computerized numerical optimization. A volume, such 
as a tumor volume (10), is graphically simulated and meshed with node 
points. The dose delivery is calculated (55, 85) depending upon input var- 
iables, deriving an objective function (56, 86) related to dose efficacy. A nu- 
merical optimization algorithm (62) optimizes the input variables based up- 
on such objective function. 
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COMPUTERIZED NUMERICAL OPTIMIZATION FOR DOSE DELIVERY 



AND TREATMENT PLANNING 



5 



CROSS-REFERENCE TO RELATED APPLICATIONS 



This application is a continuation-in-part application of U.S. 
Patent Application Serial No. 07/500,788, entitled Three -Dimensional 
Graphics Simulation and Actual Imaging Data Composite Display, to Hardy, 
filed on March 28, 1990, which is a continuation-in-part of U.S. Patent 
10 Application Serial No. 07/290,316, entitled Method and Apparatus for 

Video Presentation from a Variety of Scanner Imaging Sources, to Hardy, 
filed December 23, 1988, the teachings of both of which are incorporated 
herein by reference. 



15 No. 07/428,242, entitled Three-Dimensional Laser Localization Apparatus 
and Method for Stereotactic Diagnoses or Surgery, to Hardy, et al. , 
filed October 27, 1989, the teachings of which are incorporated herein 
by reference. 



combining three-dimensional computer graphic simulation of volumes with 
computerized numerical optimization by mathematical computation. 

Description of t he Related Art Including Information 
25 Disclosed under 37 C.F.R. fi6 1.97-l.99 ( Background Art) : 

Currently there are no efficient optimization techniques or 
treatment planning methods for optimizing the dose delivery and shape 
of, e.g., radiation zones in the treatment of various tumors in the 
body, especially the brain. Current methods involve various means of 

30 rough approximation of dose delivery to a tumor region with severe 

limitations on tailoring the radiation dose delivery to fit the exact 
dimensions and shape of the tumor volume. In addition, there is no 
method available which allows treatment planners to effectively minimize 
the dose delivery to reasonably normal brain tissue surrounding a tumor 

35 lesion and maximizing the dose to the tumor itself, other than the 



This application is also related to U.S. Patent Application Serial 



20 



BACKGROUND OF THE INVENTION 
Field of the Invention (Technical Field) : 

The present invention relates to an apparatus 



and method of 
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rather cumbersome , two-dimensional dose planning systems which rely 
heavily on dose volume histogram calculations. As designed, currently 
used systems involve very tedious iterative processing by the planner to 
fit a radiation treatment dose and volume distribution to a given tumor 
5 volume or such similar tissue volume. Current methods are very 

laborious, often requiring up to a day or more to do a single dose plan 
or any simple iterative adjustment. Such methods are severely limited 
and many dose treatment planners are constrained by current 
technological limitations to "closely approximate" a given dose 
10 delivery. There are various devices on the market; for example, the 

Leksell gamma knife (see L. Leksell, Stereotaxis and Radiosurgery: An 
operative system, Springfield, IL: Charles C. Thomas, 1971) and the 
LINAC® scalpel" (see K.R. Winston, W. Lutz, "Linear Accelerator as a 
Neurosurgical Tool for Stereotactic Radiosurgery," Neurosurgery, 
15 Vol. 22, No. 3, 1988; F. Colombo, A. Benedetti, et al. , "External 

Stereotactic Irradiation by Linear Accelerator," Neurosurgery, Vol. 16, 
No. 2, 1985; and W.A. Friedman, F.J. Bova, "The University of Florida 
Radiosurgery System," Surg. Neurol., Vol 32, pp. 334-342, 1989), which 
are capable of delivering the dose delivery to the tumor volume and 
20 optimizing the delivery such that dose zones can be tailored to tumor 

shape, configuration and size. Similarly, there are currently available 
methods for the accurate placement of radioisotopes, contained within 
catheters, in strategic areas of the brain for the direct delivery of a 
radiation dose to a given area of the brain. This method also has 
25 similar limitations in radiation dose delivery optimization (see 

K. Weaver, V. Smith, et al. , "A CT Based Computerized Treatment Planning 
System for 1-125 Stereotactic Brain Implants," Int. J', of Radiation 
Oncology, Biol. Phys . , Vol 18, pp. 445-454, 1990; P.J. Kelly, B.A. Kali, 
S. Goerss, "Computer Simulation for the Stereotactic Placement of 
30 Interstitial Radionuclide Sources into Computed Tomography -Defined Tumor 
Volumes," Neurosurgery, Vol. 14, No. 4, 1984; and B. Bauer -Kirpes , V. 
Sturm, W. Schlegel, and W.J. Lorenz, "Computerized Optimization of 1-125 
Implants in Brain Tumors," Int. J. Rjadiat. Oncol. Biol. Phys., Vol. 14 
pp. 1013-1023, 1988). 

35 None of the methods or devices described above combine three- 

dimensional graphics simulation and actual imaging data, such as 
disclosed in parent application, Serial No. 07/500,788, with 
computerized numerical optimization for dose delivery treatment 
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planning. Generally, co-pending applications S rial Nos . 07/500,788 and 
07/290,316 involve acquiring a magnetic resonance (MR), computed 
tomography (CT) , or the like, image of a volume, such as a tumor volume, 
in serial sections through the volume, outlining such sections, and 
5 creating a three-dimensional simulation of the volume. 

The concept of mathematical function optimization has existed 
since the advent of calculus and over the years has developed into very 
efficient and robust algorithms for constrained multidimensional 
nonlinear optimization. The word "optimization," when used in this 

10 context and in the specification and claim, means the rigorous use of 

algorithmic steps, implemented as computer code, to search for and find 
a mathematically defined local minimum (or maximum) of a given objective 
function. An objective function can take many forms (e.g., calculated 
stress in a structural member, aerodynamic loading on a wing, or a 

15 calculated dose of brain cell irradiation), but is simply a chosen 

measure of the desired behavior of the object, system, or process being 
designed. The term "constrained optimization" then refers to the 
optimization process, as explained above, being conducted within certain 
allowable limits or constraints. For example, a desired design 

20 objective may be to design a car frame of minimum weight design. If no 
constraints were put on this design problem, the minimum weight design 
would not be able to withstand the encountered loads during operation 
and may not even be manuf acturable . Constraints are put on the design 
problem that require the car frame to support certain loads under 

25 . various conditions and to ensure that the final design will be 

manuf acturable given current technology. Typically, "real world" design 
problems are constrained by certain necessary performance criteria. 
Modeling the physical behavior of "real world" objects, systems, and 
processes requires the use of complex nonlinear mathematical equations 

30 formed from available variables and incorporated within the computer 

code. Therefore, using the definitions provided in this paragraph, the 
term "constrained multidimensional nonlinear optimization" is defined. 

Such algorithms as discussed above were designed to replace the 
traditional "hunt and peck" process with efficient, non-random 
35 techniques for gleaning information from the computer model in the form 
of slopes and curvature of the objective function "hyper -surf ace" (a 
surface with three or more variables). When coupled in this manner, the 
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algorithms take the place of the user and autonom usly search and find 
the optimal combination of variables to maximize or minimize a desired 
ob j ective . 

Two robust and efficient nonlinear optimization algorithms 
5 available in the art are the Generalized Reduced Gradient (GRG) 

algorithm and the Sequential Quadratic Programming (SQP) algorithm. 
Since the algorithms themselves are in the public domain, many private 
versions of these algorithms are currently available as software 
packages . 

10 The present invention applies these numerical optimization 

techniques to help improve the planning process generally for optimal 
dose distribution planning and specifically for the radiation treatment 
of brain tumors. The current techniques of radiation therapy require 
the designers (neurosurgeons and radiation oncologists) to perform 

15 tedious design iterations in their search for the optimal placement of 
radioisotope seeds and external beam trajectories. By incorporating 
numerical optimization algorithms into the existing framework, the 
operation planning process becomes virtually automatic and produces 
better planning designs than current techniques and in less time. 

20 SUMMARY OF THE INVENTION 

(DISCLOSURE OF THE INVENTION') 

The preferred embodiment of the invention comprises a method and 
apparatus for optimizing dose delivery involving stereotactic computer 
techniques comprising: (a) determining an enclosing volume 

25 circumscribing a predetermined volume to be dosed; (b) meshing the 

enclosing volume with node points; (c) providing problem variables; (d) 
calculating the dose to be delivered to each of the node points; (e) 
formulating an objective function; (f) solving a numerical optimization 
algorithm minimizing the objective function; and (g) repeating (e) 

30 and (f) until the problem variables are optimized. 

The preferred method and apparatus for determining an encl sing 
volume comprises: (a) graphically simulating a calculated tumor volume 
as the predetermined volume; (b) calculating the centroid and major axis 
dimensions of the calculated tumor volume; and (c) calculating and 
35 providing an enclosing volume with the centr id of the tumor volume at 
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the center of the enclosing volume, the enclosing volume enclosing the 
tumor volume plus a nonpathologic margin. 

The preferred method and apparatus for meshing the enclosing 
volume with node points comprises: (a) assigning three-dimensional 
5 coordinates to the node points; (h) determining the position of the node 
points relative to the predetermined volume; and (c) weighting the nodes 
according to their position. 



The preferred apparatus and method for providing problem variables 
comprises (a) choosing isotopic seed activity; (b) orienting a catheter; 
(c) locating isotopic seeds along the catheter; and (d) setting upper 
and lower limits to the variables. The preferred apparatus and method 
for calculating the dose delivered to each of the node points comprises: 

(a) assigning three-dimensional coordinates to the isotopic seeds; and 

(b) calculating individual node doses by the equation: 



1 — current node number, 
j - current catheter number, 
k - current seed number, 
i2 c - total number of catheters, 
20 n s - total number of seeds per catheter, 

D ± = total radiation dose at node i, 

£ GjlSj^x^ - node radiation dose from catheter J, 

s Jk - strength of seed k in catheter j, and 

r jk ■= radius from seed k in catheter j to node i. 

25 An alternative apparatus and method for providing problem 

variables comprises : (a) choosing a pattern of beam rotation relative to 
an isocenter; (b) selecting beam strength; (c) selecting beam 
collimation; and (d) setting upper and lower limits to the problem 
variables . The preferred apparatus and method for calculating the dose 

30 delivered to each of the node points comprises: (a) separating each beam 
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path into a predetermined number of individual beams; and (b) 
calculating individual node doses by the equation: 



i — current node number, 

j - current beam number, 

n b - total number of beams , 

D ± - total radiation dose at node i, 

Gj(Sj.Cj.dj) - node radiation dose from beam j, 

s i - strength of beam J , 

Cj - collimation of beam j , and 

d~j - direction of beam j. 



The preferred apparatus and method for formulating the objective 
function comprises: (a) determining whether each node within the 
enclosing volume is dead or alive relative to a predetermined dose; 

15 (b) assigning each node to one of four specific subtotals, including (1) 
nodes inside the predetermined volume and alive, (2) nodes inside the 
predetermined volume and dead, (3) nodes outside the predetermined 
volume and alive, and (4) nodes outside the predetermined volume and 
dead; (c) dividing each of the specific subtotals by the total number of 

20 nodes; and (d) summing (1) the ratio of nodes outside the predetermined 
volume and dead and (2) the ratio of nodes inside the predetermined 
volume and alive. In this embodiment, the node error or tumor fit error 
objective function is minimized. 

An alternative apparatus and method for formulating another 
25 objective function comprises assigning status factors to the nodes as 
follows: (1) a value of +P for nodes inside the predetermined volume 
and dead, (2) a value of -I for nodes inside the predetermined volume 
and alive, (3) a value of -I for nodes outside the predetermined volume 
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and dead, and (4) a value of 0 for nodes outside the predetermined 
volume and alive, wherein: 

P - e -(Ps-W' 

Dy — Radiation dose at node i, 
D 0 - Predetermined dose, and 
5 S - decay factor. 

This embodiment for solving the numerical optimization algorithm 
comprises maximizing the node status objective function. 

A primary object of the present invention is to provide 
optimization of dose delivery in various treatments , such as treatment 
10 of brain tumors. 

A further object of the present invention is to maximize dose 
delivery in pathologic tissue while minimizing dosage in nonpathologic 
tissue. 

An advantage of the present invention is that dose delivery 
15 systems can be employed with diverse therapies. 

Yet another advantage of the present invention is the extremely 
rapid solution of the dosage delivery problem. 

Other objects, advantages and novel features, and further scope of 
applicability of the present invention will be set forth in part in the 

20 detailed description to follow, taken in conjunction with the 

accompanying drawing, and in part will become apparent to those skilled 
in the art upon examination of the following, or may be learned by 
practice of the invention. The objects and advantages of the invention 
may be realized and attained by means of the instrumentalities and 

25 combinations particularly pointed out in the appended claims. 

BRIEF DESCRIPTION OF THE DRAWINGS 
The accompanying drawings , which are incorporated into and form a 
part of the specification, illustrate several embodiments of the present 
invention and, together with the description, serve to explain the 
30 principles of the invention. The drawings are only for the purpose of 
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illustrating a preferred embodiment of the invention and are not to be 
construed as limiting the invention. 

Fig. 1 is a perspective view of a form of a dose delivery system 
showing a simulated tumor volume, brachytherapy catheters, and 
5 radioisotope seeds; 

Fig. 2 is a perspective view of a simulated enclosing volume 
circumscribing the simulated tumor volume of Fig. 1; 

Fig. 3 is a perspective view of the enclosing volume and tumor 
volume of Figs. 1 and 2 meshed with node points; 

10 Fig. 4 is a perspective view of parametric seed positioning 

relative to the enclosing volume and tumor volume of Figs. 1 and 2; 

Fig. 5 is a flow diagram of the calculation of the tumor fit error 
(objective) function of the present invention useful for brachytherapy; 

Fig. 6 is a flow chart of the optimization loop of the invention; 

15 Fig. 7 is a graph showing a curve of the "inside and dead" 

overkill penalty function of the present invention; and 

Fig. 8 is a flow diagram of the calculation of the tumor fit. error 
function of the present invention useful for external beam radiation. 

DETAILED DESCRIPTION OF T HE PREFERRED EMBODIMENTS 
20 OF THE INVENTION (BEST MODES FOR CARRYING OUT THE INVENTION) 

The present invention relates generally to the optimization of a 
dose delivery system. Specifically, the dose to be delivered can be any 
dosing agent or physical phenomena, such as radiation, poisons, chemical 
compounds, hypothermia compounds, hyperthermia- compounds , antimatter, in 

25 the form of antiprotons, photosensitive compounds, and the like. The 
invention is not limited to the dosing agents listed above, as the 
optimization techniques of the invention may be used with a variety of 
dosing agents and treatment plans. The invention is particularly useful 
for the radiation treatment of brain tumors using stereotactic 

30 techniques and thus the detailed description below is directed to such 
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treatment. Although the use of a stereotactic frame and brachy therapy- 
catheters are discussed herein, other reference frames and other dose 
delivery systems may also be used in accordance with the invention. 

Figs. 1-4 depict a tumor volume 10 calculated and graphically 
5 simulated. Methods for simulating the tumor volume 10 are disclosed in 
parent application Serial No. 07/500,788, although other simulation 
methods which might be available can also be employed. Centroid 11 and 
major axis 12 of tumor volume 10 are calculated. Depending upon the 
dose delivery system contemplated, line 13 could comprise, e.g. , the 
10 center line of brachytherapy catheter bundle 14 or the center beam of a 
major beam path employed in external beam radiation therapy. 

An enclosing volume 21 is calculated with tumor volume centroid 11 
at the center of such volume. The preferred configuration for volume 21 
is a cube, but other geometric shapes, such as a sphere, as shown in 

15 Figs. 2-4, may be employed. A volume which will enclose tumor volume 10 
is chosen such that a predetermined nonpathologic margin 32 of 
presumably healthy tissue is created. This margin is variable depending 
upon the dose delivery system employed. Normally, a higher percentage 
of healthy tissue margin is required for external beam radiation therapy 

20 to better track the dose deposited in the nonpathologic tissue. 

Subsequent to calculation of the enclosing volume, enclosing 
volume 21 is meshed with a plurality of node points, as depicted by the 
mesh intersections in Fig. 3. The number of node points is variable 
depending upon mesh dimensions, desired accuracy, and computation time 

25 available. A determination is then made as to whether each node point 
is within or without tumor volume 10. This is accomplished by 
determining the size of each voxol according to the dimension of the 
mesh and the number of nodes defined; determining the distinct 
stereotactic frame coordinates (x, j, and z) of each of these nodes; 

30 determining whether each node is located within or without tumor 
volume 10; and flagging each node accordingly (0 or I). 

Initially there are some values or variables that must be chosen 
by the user (e.g. , physician, surgeon, or technician) to start the 
optimization process. For brachytherapy or radiois topic seed 
35 implantation, the variables or values are the number of catheters to be 
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used, Che number of seeds per catheter, the configuration of the 
catheters within the available catheter block hole pattern, and the dose 
delivery time. Treatment planning may further include choosing 
radioisotope seed type {e.g., 1-125, irridium, or the like) and 
5 radioactivity level; using a nomogram to obtain a preliminary number and 
spacing of seeds relative to a given tumor volume; a computerized 
numerical optimization program to "best fit" radioisotopic seeds to the 
given tumor volume; visually inspecting graphically simulated dose zone 
delivery within the tumor volume; and confirming or adjusting the 
10 constraints of the computerized numerical optimization program. 

For external beam radiation therapy, the values or variables 
include the number of isocenters to be used and the location of the 
isocenters within the tumor. Treatment planning may further dictate 
determination of collimator size and shape and a preliminary number of 
15 isocenters for radiation beam delivery; a computerized numerical 
optimization program to "best fit" intended dose to tumor volume, 
visually inspecting graphically simulated dose zone delivery within the 
tumor volume, and confirming or adjusting the constraints of the 
computerized numerical optimization program. 

20 With specific reference to Fig. 5, the preferred initial treatment 

method for use of the invention for brachytherapy comprises the 
following steps: Tumor surface data 10 is obtained 50. This tumor 
surface data 10 is used to determine 51 an enclosing volume 21, such as 
shown in Fig. 3. The enclosing volume 21 is meshed into node points 

25 (tumor plus healthy tissue) 52, such as shown in Fig. 3. Initial 
problem variables, such as choosing the initial individual seed 
activity, the initial orientation of the catheter bundle and the initial 
location of the seeds along each catheter, are input 53. Providing the 
initial treatment plan or design may be done automatically using a 

30 nomogram incorporated into the software or input by the user. The 

initial treatment plan or design is used as a starting point for the 
optimization algorithm of the invention. The performance of the 
optimization algorithm is not all that sensitive to the initial design; 
however, it should be noted that poor initial input requires the 

35 optimization algorithm to look harder and longer for the optimal design. 
Upper and lower bounds on the design variables are either arbitrarily 
set or set to remain within available bounds (e.g. , a chosen type of 
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isotope seed will have a certain activity range in which it is 
available) . 

The isotope seed positions are then calculated and located 54 in 
x, y, and z frame coordinates. The orientation of the catheter bundle 
5 is defined by two rotation angles about the x and y axes respectively. 
Depending on the stereotactic frame design, the stereotactic frame may 
be positioned in such a way that the center of its rotation angles 
coincides or intersects with the centroid 11 of the tumor 10. This 
allows the centroid line 13 of the catheter bundle 14 to intersect with 

10 the tumor centroid 11 (see Fig. 1). An enclosing cube or other 

volume 21 is then located with its center at the tumor centroid 11 and a 
sphere with a radius just large enough to totally enclose the tumor 
volume, as in Fig. 3. Knowing the radius of the sphere and the catheter 
pattern, the intersection points of each catheter with the enclosing 

15 sphere 31 and the length of each catheter within the sphere 31 can 
easily be calculated. Since the lengths of each catheter enclosed 
within the sphere 31 will remain constant as the bundle 14 orientation 
angles are varied, these lengths can be used to parametrically locate 
the radioisotope seeds along the catheter. This is done by designating 

20 the entrance points 34 of each catheter into the enclosing volume 21 as 
a "0" in parametric space, and the exit points 36 as "I". Therefore, a 
number between 0 and I for a specific seed 38 along one of the catheters 
can be used to uniquely define its x, y, and z coordinates in the frame 
system. This concept is represented as a two-dimensional view in 

25 Fig. 4. This parametric location process using an enclosing sphere 31 
(or other volume) is implemented for two reasons: (1) to reduce the 
number of design variables (there is only one location variable for each 
seed, instead of its actual x, y, and z, coordinates), and (2) to insure 
that the relationship between the catheter bundle orientation angles and 

30 the seed positioning would remain decoupled when using a parametric 
location scheme. 



WO 91/18552 



PCT/US91/03999 



12 

Given the x, y, z coordinates of each isotope seed, the radiation 
dose to each individual node within the meshed volume can be 
calculated 55 as the sum of each seed's dose contribution to that 
particular node location (see Equation 1) . 



5 where : 

i - current node number 
J - current catheter number 
k — current seed number 
n c — total number of catheters 
10 n E - total number of seeds per catheter 

D L - total radiation dose at node i 

£ GjfSj*, Xjjj - node radiation dose from catheter J 

s^j. - strength of seed Jc in catheter J 

r Jk — radius from seed k in catheter J to node i 

15 The calculated dose received by each node in the meshed volume can be 
used in a number of ways to evaluate the "goodness" of the current 
treatment plan. 

Given the dose calculated at each node, it is determined whether 
each individual node within the meshed volume is "dead" or "alive" based 

20 on a predetermined "kill" dose of radiation. As each node 31 is 

inquired, that node is recorded as "dead" or "alive" by adding a "I" to 
the values of various subtotals. The four subtotals are identified as: 
(1) the number of nodes inside the tumor and alive; (2) the number of 
nodes inside the tumor and dead; (3) the number of nodes outside the 

25 tumor and alive; and (4) the number of nodes outside the tumor and dead. 

Using the information stored in the subtotals, the "dead" and 
"alive" ratios are calculated by dividing the number in each bin by the 
total number of nodes. These ratios are then combined to formulate 56 
an objective function for the optimization process. This objective 
30 function for tumors is referred to herein as "tumor fit error" value. 
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It is defined as th sum of the ratio f nodes outside the tumor and 
dead, plus the ratio of nodes inside the tumor and alive (see 
Equation 2) . 



Tumor Fit Error -Out d +In a 



(2) 



where : 



5 



Out d - Ratio of nodes outside the tumor and dead 
In a — Ratio of nodes inside the tumor and alive 



As shown in Fig. 6, when a value of tumor fit error is obtained 60 for 
the current seed placement and bundle orientation, a numerical 
optimization algorithm (.e.g. , Generalized Reduced Gradient or Sequential 

10 Quadratic Programming) is then used 62 to determine the optimal seed 

placement, activity and bundle orientation. Obviously, as the tumor fit 
error value approaches zero, the radiation kill zone created by the 
isotope seeds approaches the actual shape of the tumor. Therefore, the 
optimization problem can be stated: "Minimize the tumor fit error such 

15 that all nodes inside the tumor are dead." As shown in Fig. 6, the 
optimization algorithm iterates through a number of designs by 
perturbing the variable values (seed location, strength, bundle 
orientation, and other associated variables) and executing these steps 
until it reaches the best combination of variables to minimize the tumor 

20 fit error. 

With specific reference to Fig. 7, an alternative brachytherapy 
dose delivery method designed to provide homogeneous dose delivery 
comprises assigning status factors as follows: 

Calculate the individual node radiation doses in the same manner 
25 as detailed above. By comparing against the predetermined "kill" dose, 
decide if the node is "dead" or "alive" and assign it a status factor as 
defined below. 
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If the node is : 



Inside the tumor and dead, assign the node as +P 
Inside the tumor and alive, assign the node as -I 
Outside the tumor and dead, assign the node as -1 
5 Outside the tumor and alive, assign the node as 0 



P _ e ipinvf'* 

Di - Radiation dose at node i 

D 0 - Predetermined lethal radiation dose 

6 — decay factor of the overkill penalty 



10 The P function is used to define an "overkill" penalty to a node inside 
the tumor with a dose above that of the lethal dose; that is, when the 
dose level is below the lethal dose, the value of P is zero, but as it 
passes the lethal dose, its value jumps to 2 and then starts to decrease 
as the node dose level increases above the lethal limit. The 

15 parameter 5 is used to change the amount of dwell nearer to a P value 
of 1 before the exponential decay becomes effective and quickly 
decreases the "utility" of an inside node that has been overkilled. As 
can be seen, only the nodes that are inside and dead receive a factor 
greater than zero with the undesirable conditions of inside and alive 

20 and outside and dead being assigned negative status factors. 

After assigning a status factor to each node, an objective 
function for the optimization process can be formed by summing all of 
the node status factors (see Equation 3). This "tumor status" objective 
function will be different for each seed placement design and would 

25 reach its maximum when there are few "outside dead" and "inside alive" 
nodes and most of the "inside dead" nodes have a status factor of I 
(meaning they have just attained a level of "dead"). Notice that 
maximizing the "tumor status" function will have the effect of 
homogenizing the dose distribution within the tumor boundaries because 

30 the maximum value of the node status factor is assigned only to those 
nodes that are not "overkilled" by too much. Therefore, the 
optimization problem can be stated: "Maximize the tumor status 
function." Insuring that all the nodes inside are dead could be 
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accomplished by increasing th negative value of the "inside alive" 
status factor and more heavily penalizing the "tumor status" function 
when -inside alive" nodes exist. The optimization algorithm will 
iterate through a number of designs or treatment plans by perturbing the 
5 variable values (seed location, strength, bundle orientation, and other 
associated variables) until it reaches the best combination of variables 
to maximize the "tumor status" function. 

Employment of external beam radiotherapy also requires an initial 
design or treatment plan comprising problem variables. With specific 

10 reference to Fig. 8, and as in Fig. 5, tumor surface data is 

obtained 80, which is used to determine 81 an enclosing volume. Initial 
program variables are chosen and input 82. These program variables 
include a pattern of beam rotation for each isocenter, the beam strength 
and beam collimation (either or both of the last two could be provided 

15 as a function of beam rotation with the variables then becoming the 
coefficients of polynomial functions or control points of Bezier 
curves). An initial treatment plan appropriate for this type of therapy 
is used as a starting point for the optimization algorithm. The 
performance of the optimization algorithm is not all that sensitive to 

20 the initial treatment plan; however, it should be noted that poor 

initial input will again require the optimization algorithm to look 
harder and longer for the optimal design. Upper and lower bounds on the 
design variables are either arbitrarily set or set to remain within 
available bounds (e.g. , a beam rotation pattern may be constrained. to 

25 remain between certain angles to avoid sensitive areas of the brain) . 

Each major beam path within each isocenter 's application pattern 
is separated 84 into a predetermined number of individual beams. This 
is done so that the actual continuous beam paths can be approximated and 
evaluated by the dose function. 

30 Given the orientation angles (direction cosines) of each discrete 

beam, the radiation dose to each individual node is calculated 85 within 
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the meshed volume by summing each beam's dose contribution to each node 
as follows : 

where : 

i — current node number 
5 J — current beam number 

n b — total number of beams 

D ± - total radiation dose at node i 

G^Sj. c j . dj) - node radiation dose from beam J 

Sj — strength of beam J 
10 Cj - collimation of beam j 

d~j - direction of beam J 

The objective function {e.g., tumor fit function) is formulated 86 
as discussed above. Given the dose calculated at each node, it is 
determined whether each individual node within the meshed volume is 

15 "dead" or "alive" based on a predetermined "kill" dose of radiation. As 
each node is inquired, that node is recorded as "dead" or "alive" by 
adding a "1" to the values of various collection subtotals. The four 
subtotals identified are: (1) the number of nodes inside the tumor and 
alive; (2) the number of nodes inside the tumor and dead; (3) the number 

20 of nodes outside the tumor and alive; and (4) the number of nodes 
outside the tumor and dead. 

Using the information stored in the subtotals, the "dead" and 
"alive" ratios are calculated by dividing the number in each subtotal by 
the total number of nodes. These ratios are then combined to form an 
25 objective function for the optimization process. This objective 

function for tumors is referred to herein as the "tumor fit error" 
value. It is defined as the sum of the ratio of nodes outside the tumor 
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and dead, plus the ratio of nodes insid the tumor and alive (see 
Equation 5) . 

Tumor Fi t Error -Ou t d + In B ( 5 ) 



Ou£ d - Ratio of nodes outside the tumor and dead 
Jn a - Ratio of nodes inside the tumor and alive 

When a tumor fit error value is obtained for the given beam variables 
and isocenters, a numerical optimization algorithm (e.g., Generalized 
Reduced Gradient or Sequential Quadratic Programming algorithms) is used 
to determine the optimal combination of beam variables. As the tumor 
fit error value approaches zero the radiation kill zone created by the 
beam trajectories approaches the actual shape of the tumor. Therefore, 
the optimization problem can be stated: "Minimize the tumor fit error 
such that all nodes inside the tumor are dead." Again, as shown in 
Fig. 6, the optimization algorithm iterates through a number of designs 
by perturbing the variable values (beam rotation pattern, strength, and 
collimation, and other associated variables) until it reaches the best 
combination of variables to minimize the tumor fit error. 



Alternatively, a homogeneous dose delivery method can also be used 
with external beam radiotherapy using status factors, as previously 
20 discussed. 

As discussed above, the computerized numerical optimization dose 
delivery methods herein disclosed are also suitable for other dose 
delivery methods and devices including thermal or cryogenic probes, 
infrared or laser systems, chemical doses, drug diffusion systems, and, 
25 ultimately, antimatter in the form of antiprotons . 

Additionally, other modifications may be made to the optimal dose 
delivery system of the invention. These are discussed below. 

Weighting of nodes in sensitive brain areas may exert greater 
influence on the objective function, ther by lessening dosage on such 
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nodes . Such capability could be directly tied to brain nap databases , 
such as disci sed in co-pending application Serial No. 07/290,316. 

In an alternative embodiment, increasing node density near the 
tumor boundaries while decreasing node density near the tumor center and 
5 boundaries of the enclosing volume also improves "best fit" of kill zone 
to tumor volume, thereby increasing the accuracy of the treatment. 

Use of artificial intelligence rule bases in conjunction with the 
optimization algorithms allows more flexibility in defining the 
optimization process. For example, artificial intelligence can be used 

10 to allow the number of catheters and their layout pattern to be 

variables in brachy therapy optimization. This occurs because the actual 
solution method changes as the number of catheters and their pattern 
changes. This integration of artificial intelligence into the numerical 
optimization process allows current parameters which are entered as 

15 constants to the equations to be given a desired range ("limited 
variables") . 

Real-time three dimensional (3-D) visualization of the 
optimization process can be incorporated. The radiation treatment 
planner can watch the dose distribution change to fit the tumor volume. 
20 This permits even greater flexibility in the planning process, allowing 
the treatment planner to stop the process at any point during the 
optimization, possibly adjust "problem variables" as appropriate, and 
re-start the optimization process with the new entries. 

Three dimensional simulation of dose-volume histograms can be 
25 presented to the user in a fashion such that the user can simultaneously 
perceive the tumor volume in relationship to a 3-D rendition of a dose- 
volume histogram of dosimetry to surrounding tissue. 

Currently there are no highly accurate methods of determining 
actual tumor boundaries ; use of radioisotope or other such tagging or 
30 identification techniques in conjunction with imaging and graphics 
techniques, such as disclosed in co-pending applications, helps to 
better define the actual tumor margin. 
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Although the invention has be n described with reference to these 
preferred embodiments, other embodiments can achieve the same results. 
Variations and modifications of the present invention will be obvious to 
those skilled in the art and it is intended to cover in the appended 
5 claims all such modifications and equivalents. 
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What is claimed is: 

1. A method of optimizing dose delivery involving stereotactic 
computer techniques comprising the steps of: 

a. determining an enclosing volume circumscribing a 
predetermined volume to be dosed; 

meshing said enclosing volume with node points; 
providing problem variables; 

calculating the dose to be delivered to each of said 



node points; 



the objective 
are optimized 



formulating an objective function; 

solving a numerical optimization algorithm minimizing 
function; and 

repeating steps e. and f. until said problem variables 



2. The method of claim 1 wherein the step of determining an 
enclosing volume circumscribing the predetermined volume to be dosed 
comprises the steps of: 

a. graphically simulating a calculated tumor volume as 
said predetermined volume; 

b. calculating the centroid and major axis dimension of 
said calculated tumor volume; and 

c. calculating and providing an enclosing volume with 
said centroid of said tumor volume at the center of the enclosing 
volume, said enclosing volume enclosing said tumor volume plus a 
nonpathologic margin. 

3. The method of claim 1 wherein the step of meshing said 
enclosing volume with node points comprises the steps of: 

a. assigning three-dimensional coordinates to said node 

points ; 

b. determining the position of said node points relative 
to said predetermined volume; and 

c. weighting said nodes according to their position. 
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4. The method f claim 1 wherein said step of providing problem 
variables comprises the steps of: 

a. choosing isotopic seed activity; 

b. orienting a catheter; 

5 c. locating isotopic seeds along the catheter; and 

d. setting upper and lower limits to said variables. 



5. The method of claim 4 wherein said step of calculating the 
dose delivered to each of said node points comprises the steps of: 

a. assigning three-dimensional coordinates to the 
isotopic seeds; and 

b . calculating individual node doses by the equation 



1 - current node number, 
j — current catheter number, 
15 k - current seed number, 

n c — total number of catheters, 

n B — total number of seeds per catheter, 

D L — total radiation dose at node i , 

52 Gj(S Jfc ,r Jfc ) - node radiation dose from catheter j, 

20 s jk - strength of seed k in catheter j , and 

r Jk - radius from seed k in catheter J to node i. 

6 . The method of claim 1 wherein said step of providing problem 
variables comprises the steps of: 

a. choosing a pattern of beam rotation relative to an 

25 isocenter; 

b. selecting beam strength; 

c. selecting beam collimation; and 

d. setting upper and lower limits to said problem 



variables . 
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7. The method of claim 6 wherein Che step of calculating the 
dose delivered to each of said node points comprises the steps of: 

a. separating each beam path into a predetermined number 
of individual beams; and 

b. calculating individual node doses by the equation 



i - current node number, 

J - current beam number, 

n b - total number of beams , 

Dj_ - total radiation dose at node i, 

G^Sj.Cj.dj) — node radiation dose from beam J, 

Sj - strength of beam j, 

Cj - collimation of beam j , and 

d } - direction of beam J . 



15 8. The method of claim 1 wherein the step of formulating an 

objective function comprises the steps of: 

a. determining whether each said node within said 
enclosing volume is dead or alive relative to a predetermined dose; 

b. assigning each said node to one of four specific 
20 subtotals, including (1) nodes inside the predetermined volume and 

alive, (2) nodes inside the predetermined volume and dead, (3) nodes 
outside the predetermined volume and alive, and (4) nodes outside the 
predetermined volume and dead; 

c. dividing each of said specific subtotals by the total 
25 number of nodes ; and 

d. summing (1) the ratio of nodes outside the 
predetermined volume and dead and (2) the ratio of nodes inside the 
predetermined volume and alive. 

9 . The method of claim 8 wherein the step of formulating an 
30 objective function comprises minimizing the objective function. 
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10. The method of claim 1 wherein the step of formulating an 
objective function comprises assigning status factors to the nodes as 
follows: (1) a value of +P for nodes inside the predetermined volume 
and dead, (2) a value of -2 for nodes inside the predetermined volume 

5 and alive, (3) a value of -1 for nodes outside the predetermined volume 
and dead, and (4) a value of 0 for nodes outside the predetermined 
volume and alive, wherein: 

p _ e HPi-W 

D y - Radiation dose at node i, 
D 0 - Predetermined dose, and 
10 6 - decay factor. 

11 . The method of claim 10 wherein the step of solving a 
numerical optimization algorithm comprises maximizing the objective 
function. 

12. A apparatus for optimizing dose delivery involving 
stereotactic computer techniques comprising: 

means for determining an enclosing volume circumscribing a 
predetermined volume to be dosed; 

means for meshing said enclosing volume with node points; 
means for providing problem variables ; 

means for calculating the dose to be delivered to each, of 
said node points; 

means for formulating an objective function; and 
means for solving a numerical optimization algorithm 
minimizing the objective function. 



15 



20 
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13. The apparatus of claim 12 wherein said means for determining 
an enclosing volume circumscribing Che predetermined v lume to be dosed 
comprises : 

means for graphically simulating a calculated tumor volume 
as said predetermined volume; 

means for calculating the centroid and major axis dimension 
of said calculated tumor volume; and 

means for calculating and providing an enclosing volume with 
said centroid of said tumor volume at the center of the enclosing 
volume, said enclosing volume enclosing said tumor volume plus a 
nonpathologic margin. 

14. The apparatus of claim 12 wherein said means for meshing 
said enclosing volume with node points comprises: 

means for assigning three-dimensional coordinates to said 

node points; 

means for determining the position of said node points 
relative to said predetermined volume; and 

means for weighting said nodes according to their position. 

15 . The apparatus of claim 12 wherein said means for providing 
problem variables comprises: 

means for choosing isotopic seed activity; 
means for orienting a catheter; 

means for locating isotopic seeds along the catheter; and 
means for setting upper and lower limits to said variables. 
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16. The apparatus of claim 15 wherein said means for calculating 
the dose delivered to each of said node points comprises : 

means for assigning three-dimensional coordinates to the 
isotopic seeds; and 

means for calculating individual node doses by the equation 

wherein: 

i - current node number, 
j - current catheter number, 
k — current seed number, 
10 n c - total number of catheters, 

n s - total number of seeds per catheter, 
D t - total radiation dose at node i, 

Gj(S Jk ,z lk ) - node radiation dose from catheter J, 

s Jk — strength of seed k in catheter J, and 
15 r Jk - radius from seed k in catheter J to node i. 

17 . The apparatus of claim 12 wherein said means for providing 
problem variables comprises: 

means for choosing a pattern of beam rotation relative to an 

isocenter; 

20 means for selecting beam strength; 

means for selecting beam collimation; and 

means for setting upper and lower limits to said problem 

variables . 
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18. The apparatus of claim 17 wherein said means for calculating 
the dose delivered to each of said node points comprises: 

means for separating each beam path into a predetermined 
number of individual beams; and 
5 means for calculating individual node doses by the equation 



i - current node number, 

J — current beam number , 

n b - total number of beams, 

B l — total radiation dose at node i , 

G^Sj. Cj, dj) - node radiation dose from beam j, 

Sj - strength of beam j , 

c i - collimation of beam J, and 

d~j — direction of beam J. 



19. The apparatus of claim 12 wherein said means for formulating 
an objective function comprises: 

means for determining whether each said node within said 
enclosing volume is dead or alive relative to a predetermined dose; 

means for assigning each said node to one of four specific 
subtotals, including (1) nodes inside the predetermined volume and 
alive, (2) nodes inside the predetermined volume and dead, (3) nodes 
outside the predetermined volume and alive, and (4) nodes outside the 
predetermined volume and dead; 

means for dividing each of said specific subtotals by the 
total number of nodes; and 

means for summing (1) the ratio of nodes outside the 
predetermined volume and dead and (2) the ratio of nodes inside the 
predetermined volume and alive. 



20. The method of claim 19 wherein the step of formulating an 
30 objective function comprises minimizing the objective function. 
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21. The apparatus of claim 12 wherein said means for formulating 
an objective function comprises assigning status factors to the nodes as 
follows: (1) a value of +P for nodes inside the predetermined volume 
and dead, (2) a value of -2 for nodes inside the predetermined volume 

5 and alive, (3) a value of -2 for nodes outside the predetermined volume 
and dead, and (4) a value of 0 for nodes outside the predetermined 
volume and alive, wherein: 

p _ e -<Pi-W 

- Radiation dose at node i, 
D 0 — Predetermined dose, and 
10 6 — decay factor. 

22. The apparatus of claim 21 wherein said means for solving a 
numerical optimization algorithm comprises maximizing the objective 
function. 
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